RESULTS excluding rare and separating species by groups
Running models with different datasets:
all: original data, all species all individuals
exclude: exclude rare species
regroup: regrouping observations of rare species in ONE species group
Using only the 5x5 m quadrat scale
#Figures colors/labels
cores <- c( quadrat = "#c00000",
`sp:time`= "#203864",
time = "#ed7d31",
`quadrat:time` = "#843c0c",
`quadrat:sp`= "#b666d2",
sp = "#70ad47",
Residual="#ffc000")
labvit <- as_labeller(c(grow = "Growth", mort="Mortality",
rec="Recruitment"))
grpTlab = as_labeller(c(quadrat = "space", `sp:time`="species x time",
`quadrat:time`="space x time", `quadrat:sp`= "species x space",
sp="species", Residual = "residual",
time= "time"))
grpvitT <- as_labeller(c(grow = "Growth", mort="Mortality",
rec="Recruitment",
quadrat = "space", `sp:time`="species*time",
`quadrat:time`="space*time", `quadrat:sp`= "species*space",
sp="species", Residual = "residual",
time= "time"))
labplot = as_labeller(c(quadrat = "space",
`sp:time`="species x time",
`quadrat:time`="space x time",
`quadrat:sp`= "species x space",
sp="species",
Residual = "residual",
time= "time",
fus="Fushan", luq = "Luquillo", lam="Lambir",
pas = "Pasoh", bci="Barro Colorado"))Data
Growth
grp <- c("all","exclude","regroup")
path1 = here("models_outputs/models_time/")
tudo <- list()
for (i in 1:length(grp)){
path <- paste0(path1, grp[i], "/grow/table")
files = list.files(path)
# excluding wyw
files = files[grep("wyw", files, invert = T)]
all <- lapply(files, function(x) get(assign(x, load(paste0(path, "/", x)))))
names(all) <- substr(files,1,nchar(files)-6)
tudo[[i]] <- all
}
names(tudo) <- grp
tes <- map(tudo,bind_rows, .id="id") %>% bind_rows() %>% filter(q_size=="quad_5")
#tes %>% count(data,fplot) %>% mutate(n=n/8)
table(tes$data,tes$fplot)/8##
## bci fus lam luq pas
## all 10 10 10 10 10
## exclude 10 10 10 10 10
## regroup 10 10 10 10 10
#remove intercept to include again as a separate column
intercepts <- tes %>% filter(term == "intercept")
names(intercepts)[9] = "intercept"
grow <- tes %>% filter(term != "intercept") %>%
left_join(intercepts[,c(1:4,9)], by= c("id", "sub","data", "fplot")) %>%
group_by(id, sub, data, fplot) %>%
mutate(VPC = variance/sum(variance)) %>% ungroup() %>%
mutate(term = fct_relevel(term,"time", "quadrat", "sp:time", "quadrat:time",
"quadrat:sp", "sp", "Residual")) %>%
as.data.frame()
# divide the sd of the terms by the mean growth (intercept) for each dataset.
grow$stdev = grow$Estimate/grow$interceptMort
tudo <- list()
for (i in 1:length(grp)){
path <- paste0(path1,grp[i], "/mort/table")
files = list.files(path) # excluding wyw
files = files[grep("wyw", files, invert = T)]
all <- lapply(files, function(x) get(assign(x, load(paste0(path, "/", x)))))
names(all) <- substr(files,1,nchar(files)-6)
tudo[[i]] <- all
}
names(tudo) <- grp
tes <- map(tudo,bind_rows, .id="id") %>% bind_rows() %>% filter(q_size=="quad_5")
#tes %>% count(data,fplot) %>% mutate(n=n/8)
table(tes$data,tes$fplot)/8##
## bci fus lam luq pas
## all 10 10 10 10 10
## exclude 10 10 10 10 10
## regroup 10 10 10 10 10
#remove intercept to include again as a separate column
intercepts <- tes %>% filter(term == "intercept")
names(intercepts)[9] = "intercept"
mort <- tes %>% filter(term != "intercept") %>%
left_join(intercepts[,c(1:4,9)], by= c("id", "sub","data", "fplot")) %>%
group_by(id, sub, data, fplot) %>%
mutate(VPC = variance/sum(variance)) %>% ungroup() %>%
mutate(term = fct_relevel(term,"time", "quadrat", "sp:time", "quadrat:time",
"quadrat:sp", "sp", "Residual")) %>%
as.data.frame()
mort$stdev <- mort$EstimateRec
tudo <- list()
for (i in 1:length(grp)){
path <- paste0(path1, grp[i], "/rec/table")
files = list.files(path) # excluding wyw
files = files[grep("wyw", files, invert = T)]
all <- lapply(files, function(x) get(assign(x, load(paste0(path, "/", x)))))
names(all) <- substr(files,1,nchar(files)-6)
tudo[[i]] <- all
}
names(tudo) <- grp
tes <- map(tudo,bind_rows, .id="id") %>% bind_rows() %>% filter(q_size=="quad_5")
#tes %>% count(data,fplot) %>% mutate(n=n/8)
table(tes$data,tes$fplot)/8##
## bci fus lam luq pas
## all 10 10 10 10 10
## exclude 10 10 10 10 10
## regroup 10 10 10 10 10
#remove intercept to include again as a separate column
intercepts <- tes %>% filter(term == "intercept")
names(intercepts)[9] = "intercept"
rec <- tes %>% filter(term != "intercept") %>%
left_join(intercepts[,c(1:4,9)], by= c("id", "sub","data", "fplot")) %>%
group_by(id, sub, data, fplot) %>%
mutate(VPC = variance/sum(variance)) %>% ungroup() %>%
mutate(term = fct_relevel(term,"time", "quadrat", "sp:time", "quadrat:time",
"quadrat:sp", "sp", "Residual")) %>%
as.data.frame()
rec$stdev <- rec$EstimateCombining all results and taking the mean stdev and VPC among subsamples
NOT taking into account the same dataset (sub)
comb <- bind_rows(list(grow=grow, mort=mort, rec=rec), .id="vital") %>%
mutate(term = fct_relevel(term, "quadrat", "sp:time","time", "quadrat:time",
"quadrat:sp", "sp","Residual")) %>%
ungroup() %>%
mutate(data = fct_relevel(data, "all","regroup", "exclude"))
mcomb <- comb %>% group_by(vital,data,fplot,q_size,term) %>%
summarise(stdev = mean(stdev),
VPC = mean(VPC),
richness = round(mean(richness))) %>%
filter(q_size == "quad_5")number of fplots with results
table(mcomb$vital, mcomb$data)/7##
## all regroup exclude
## grow 5 5 5
## mort 5 5 5
## rec 5 5 5
INCLUDING THE % of rare species
load(here("data/rare_species.Rdata"))
comb <- comb %>% left_join(rare, by="fplot")
mcomb <- mcomb %>% left_join(rare, by="fplot")Comparing groups
mcomb %>%
ggplot(aes(x=data, y=VPC)) +geom_col(aes(fill=term, color=term), alpha=0.7) +
scale_fill_manual(values=cores) +
scale_color_manual(values=cores) +
facet_grid(vital~fplot) +
theme(axis.text.x = element_text(angle=45, hjust=1))fig <- comb %>% filter(vital=="grow") %>%
ggplot(aes(x=data, y=VPC, col=term)) +
geom_quasirandom(alpha=0.5, groupOnX = TRUE) +
facet_grid(fplot~term, labeller=labplot) +
stat_summary(fun=median, size=0.3, fatten=1,col="black",
geom="crossbar", width=0.7)+
scale_color_manual(values = cores) +
xlab("") +
theme(legend.position = "none",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45,hjust=1)) +
labs(title= "VPC per subsample of 5ha - GROWTH")
figfig <- comb %>% filter(vital=="mort") %>%
ggplot(aes(x=data, y=VPC, col=term)) +
geom_quasirandom(alpha=0.5, groupOnX = TRUE) +
facet_grid(fplot~term, labeller=labplot) +
stat_summary(fun=median, size=0.3, fatten=1,col="black",
geom="crossbar", width=0.7)+
scale_color_manual(values = cores) +
xlab("") +
theme(legend.position = "none",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45,hjust=1)) +
labs(title= "VPC per subsample of 5ha - MORT")
figfig <- comb %>% filter(vital=="rec") %>%
ggplot(aes(x=data, y=VPC, col=term)) +
geom_quasirandom(alpha=0.5, groupOnX = TRUE) +
facet_grid(fplot~term, labeller=labplot) +
stat_summary(fun=median, size=0.3, fatten=1,col="black",
geom="crossbar", width=0.7)+
scale_color_manual(values = cores) +
xlab("") +
theme(legend.position = "none",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45,hjust=1)) +
labs(title= "VPC per subsample of 5ha - REC")
figmean over fplots
mcomb %>% group_by(vital,data,term) %>% summarise(VPC=mean(VPC)) %>%
ggplot(aes(x=data, y=VPC)) +geom_col(aes(fill=term, color=term), alpha=0.7) +
scale_fill_manual(values=cores) +
scale_color_manual(values=cores) +
facet_grid(~vital)+
theme(axis.text.x = element_text(angle=45, hjust=1))comb %>%
ggplot(aes(x=data,y=VPC, col=data)) + geom_boxplot()+
facet_grid(vital~term) +
theme(axis.text.x = element_text(angle=45,hjust=1),
legend.position = "none")mcomb %>% mutate(xis = as.numeric(data)) %>%
ggplot(aes(x=xis,y=stdev, col=fplot)) + geom_line() +
facet_grid(vital~term, scales = "free") +
scale_x_continuous(breaks = 1:3, labels=levels(mcomb$data),
name="") +
theme(axis.text.x = element_text(angle=45, hjust=1))USING the same dataset for the comparisons
It means, for each forest, linking the same subsampled data used
fig <- comb %>% filter(vital=="grow") %>%
ggplot(aes(x=as.numeric(as.factor(data)), y=VPC, col=sub)) +
geom_line() +
facet_grid(fplot~term, labeller=labplot) +
scale_x_continuous(breaks=1:3, labels=levels(comb$data))+
xlab("") +
theme(legend.position = "none",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45,hjust=1)) +
labs(title= "VPC per subsample of 5ha - GROW")
figfig <- comb %>% filter(vital=="mort") %>%
ggplot(aes(x=as.numeric(as.factor(data)), y=VPC, col=sub)) +
geom_line() +
facet_grid(fplot~term, labeller=labplot) +
scale_x_continuous(breaks=1:3, labels=levels(comb$data))+
xlab("") +
theme(legend.position = "none",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45,hjust=1)) +
labs(title= "VPC per subsample of 5ha - mort")
figfig <- comb %>% filter(vital=="rec") %>%
ggplot(aes(x=as.numeric(as.factor(data)), y=VPC, col=sub)) +
geom_line() +
facet_grid(fplot~term, labeller=labplot) +
scale_x_continuous(breaks=1:3, labels=levels(comb$data))+
xlab("") +
theme(legend.position = "none",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45,hjust=1)) +
labs(title= "VPC per subsample of 5ha - rec")
figlabplot2 = as_labeller(c(grow = "Growth", mort="Mortality",
rec="Recruitment",
fus="Fushan", luq = "Luquillo", lam="Lambir",
pas = "Pasoh", bci="Barro Colorado"))
fg <- mcomb %>% filter(data %in% c("all", "exclude", "regroup")) %>%
mutate(data = fct_relevel(data, "all", "exclude", "regroup")) %>%
ggplot(aes(x=data, y=VPC)) +
geom_col(aes(fill=term, color=term), alpha=0.7) +
scale_fill_manual(values=cores,labels=
c("Space","Species x time", "Time", "Space x time", "Species x Space","Species","Residual")) +
scale_color_manual(values=cores, labels=
c("Space","Species x time", "Time", "Space x time", "Species x Space","Species","Residual")) +
facet_grid(fplot~vital, labeller = labplot2) +
theme(axis.text.x = element_text(angle=45, hjust=1))
fg# png("figs/exclude_regroup_rare_VPC_time_bar2.png", height = 800, width=800, res=100)
# fg
# dev.off()mcomb %>% filter(data %in% c("all", "exclude", "regroup")) %>%
mutate(x=as.factor(data)) %>% mutate(x= as.numeric(x)) %>%
ggplot(aes(x=x, y=stdev, col=fplot)) + geom_line() +
facet_grid(vital~term, scales="free") +
scale_x_continuous(breaks = 1:3, labels=c("all","regroup", "exclude"),
name="")+
theme(axis.text.x = element_text(angle=45, hjust=1),
legend.position = "right") +
ggtitle("SD")mcomb %>% filter(data %in% c("all","regroup", "exclude")) %>%
mutate(x=as.factor(data)) %>% mutate(x= as.numeric(x)) %>%
ggplot(aes(x=data, y=stdev, col=data)) + geom_boxplot() +
facet_grid(vital~term, scales="free") +
theme(axis.text.x = element_text(angle=45, hjust=1),
legend.position = "none") +
ggtitle("SD")mcomb %>% filter(data %in% c("all","regroup", "exclude")) %>%
mutate(x=as.factor(data)) %>% mutate(x= as.numeric(x)) %>%
ggplot(aes(x=x, y=VPC, col=fplot)) + geom_line() +
facet_grid(vital~term, scales="free") +
scale_x_continuous(breaks = 1:3, labels=c("all","regroup", "exclude"),
name="")+
theme(axis.text.x = element_text(angle=45, hjust=1),
legend.position = "right") +
ggtitle("VPC")mcomb %>% filter(data %in% c("all","regroup", "exclude")) %>%
mutate(x=as.factor(data)) %>% mutate(x= as.numeric(x)) %>%
ggplot(aes(x=data, y=VPC, col=data)) + geom_boxplot() +
facet_grid(vital~term, scales="free") +
theme(axis.text.x = element_text(angle=45, hjust=1),
legend.position = "none") +
ggtitle("VPC")Difference to all data
Absolute difference
dif <- mcomb %>% filter(data %in% c("all","regroup", "exclude")) %>%
mutate(data = fct_relevel(data, "all", "exclude", "regroup")) %>%
pivot_wider(id_cols = c(vital,fplot,term, p.rare),
names_from = data,
values_from = c(stdev,VPC, richness)) %>%
mutate(stdev_difr.regroup = stdev_regroup - stdev_all,
stdev_difr.exclude = stdev_exclude - stdev_all,
VPC_difr.regroup = VPC_regroup - VPC_all,
VPC_difr.exclude = VPC_exclude - VPC_all,
prop.rare.sp = (richness_exclude/richness_all)) %>%
select(1:4,14:18) %>%
pivot_longer(5:8, names_to = c(".value", "group"),
names_pattern = "(.*)_(.*)") %>%
mutate(group = substr(group,5,nchar(group)))dif %>%
ggplot(aes(x=group, y=stdev)) +
facet_grid(vital~term, scales= "free_x") +
geom_boxplot() +
geom_quasirandom(aes(col=fplot)) +
ylab("SD group - SD all") + xlab("") +
geom_hline(yintercept = 0, linetype= "dotted") +
ggtitle("Absolute difference in SD") +
theme(legend.position = "bottom",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45, hjust=1))dif %>%
ggplot(aes(x=group, y=VPC)) +
facet_grid(vital~term, scales= "free_x") +
geom_boxplot(outlier.alpha = 0) +
geom_quasirandom(aes(col=fplot)) +
ylab("VPC group - VPC all") + xlab("") +
geom_hline(yintercept = 0, linetype= "dotted") +
ggtitle("Absolute difference in VPC") +
theme(legend.position = "bottom",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45, hjust=1))resdif <- dif %>% group_by(vital, term, group) %>% summarise(mean=mean(VPC, na.rm=T)) %>%
pivot_wider(names_from = term, values_from = mean)Relative differences in SD to all data
(group value - all data value)*100 / all data value
difr <- mcomb %>% filter(data %in% c("all","regroup", "exclude")) %>%
mutate(data = fct_relevel(data, "all", "exclude", "regroup")) %>%
pivot_wider(id_cols = c(vital,fplot,term, p.rare),
names_from = data,
values_from = c(stdev,VPC, richness)) %>%
mutate(stdev_difr.regroup = (stdev_regroup - stdev_all)*100/stdev_all,
stdev_difr.exclude = (stdev_exclude - stdev_all)*100/stdev_all,
VPC_difr.regroup = (VPC_regroup - VPC_all)*100/VPC_all,
VPC_difr.exclude = (VPC_exclude - VPC_all)*100/VPC_all,
prop.rare.sp = (richness_exclude/richness_all)) %>%
select(1:4,14:18) %>%
pivot_longer(5:8, names_to = c(".value", "group"),
names_pattern = "(.*)_(.*)") %>%
mutate(group = substr(group,5,nchar(group)))difr %>%
ggplot(aes(x=group, y=stdev)) +
facet_grid(vital~term, scales= "free_x") +
geom_boxplot() +
geom_quasirandom(aes(col=fplot)) +
ylab("SD proportional diference (%)") + xlab("") +
geom_hline(yintercept = 0, linetype= "dotted") +
ggtitle("Relative difference in SD") +
theme(legend.position = "bottom",
panel.background = element_rect(colour="lightgray"),
axis.text.x = element_text(angle=45, hjust=1))resdifr <- difr %>% group_by(vital, term, group) %>% summarise(mean=mean(stdev, na.rm=T)) %>%
pivot_wider(names_from = term, values_from = mean)
resdifr %>% kable(digits = 2)| vital | group | quadrat | sp:time | time | quadrat:time | quadrat:sp | sp | Residual |
|---|---|---|---|---|---|---|---|---|
| grow | .exclude | 11.65 | -13.84 | 3.55 | -3.95 | 1.86 | -11.99 | 1.53 |
| grow | .regroup | 13.69 | -23.98 | 7.45 | -6.76 | 0.31 | -12.50 | 8.09 |
| mort | .exclude | 5.47 | -3.67 | -0.48 | 0.30 | -3.07 | -5.19 | 0.00 |
| mort | .regroup | 6.95 | -5.41 | 0.07 | -2.83 | -0.94 | -6.54 | 0.00 |
| rec | .exclude | 21.03 | -7.91 | 6.51 | 1.16 | 2.38 | -7.74 | 0.00 |
| rec | .regroup | 24.34 | -10.11 | 5.37 | 1.79 | -3.80 | -10.61 | 0.00 |
Relative differences in SD
difr %>% group_by(vital, term, group) %>% summarise(mean=mean(stdev, na.rm=T)) %>%
pivot_wider(names_from = c(vital,group), values_from = mean) %>%
kable(digits=1)| term | grow_.exclude | grow_.regroup | mort_.exclude | mort_.regroup | rec_.exclude | rec_.regroup |
|---|---|---|---|---|---|---|
| quadrat | 11.7 | 13.7 | 5.5 | 7.0 | 21.0 | 24.3 |
| sp:time | -13.8 | -24.0 | -3.7 | -5.4 | -7.9 | -10.1 |
| time | 3.5 | 7.5 | -0.5 | 0.1 | 6.5 | 5.4 |
| quadrat:time | -3.9 | -6.8 | 0.3 | -2.8 | 1.2 | 1.8 |
| quadrat:sp | 1.9 | 0.3 | -3.1 | -0.9 | 2.4 | -3.8 |
| sp | -12.0 | -12.5 | -5.2 | -6.5 | -7.7 | -10.6 |
| Residual | 1.5 | 8.1 | 0.0 | 0.0 | 0.0 | 0.0 |
Absolute differences in VPC
dif %>% group_by(vital, term, group) %>% summarise(mean=mean(VPC, na.rm=T)) %>%
pivot_wider(names_from = c(vital,group), values_from = mean) %>%
kable(digits=3)| term | grow_.exclude | grow_.regroup | mort_.exclude | mort_.regroup | rec_.exclude | rec_.regroup |
|---|---|---|---|---|---|---|
| quadrat | 0.002 | 0.002 | 0.000 | 0.002 | 0.001 | 0.002 |
| sp:time | -0.007 | -0.008 | -0.001 | -0.001 | -0.005 | -0.006 |
| time | 0.002 | 0.002 | -0.001 | 0.000 | 0.009 | 0.011 |
| quadrat:time | 0.005 | 0.002 | 0.005 | 0.000 | 0.003 | 0.013 |
| quadrat:sp | 0.014 | -0.006 | -0.002 | -0.001 | 0.007 | -0.009 |
| sp | -0.042 | -0.050 | -0.020 | -0.024 | -0.021 | -0.026 |
| Residual | 0.026 | 0.059 | 0.019 | 0.025 | 0.006 | 0.016 |
Comparing with Proportion of rare species
dif %>%
ggplot(aes(x=p.rare, y=stdev, linetype=group)) +
geom_point(aes(col=fplot)) +
geom_smooth(method="lm", se=F) +
facet_grid(vital~term) +
ylab("Difference in SD to original data") +
xlab("Prop of rare species") +
geom_hline(yintercept = 0, linetype= "dotted") +
ggtitle("Difference in SD") +
theme(axis.text.x = element_text(angle=45, hjust=1))difr %>%
ggplot(aes(x=p.rare, y=stdev, linetype=group)) +
geom_point(aes(col=fplot)) +
geom_smooth(method="lm", se=F) +
facet_grid(vital~term) +
ylab("Proportional difference in SD to original data") +
xlab("Proportion of rare species") +
geom_hline(yintercept = 0, linetype= "dotted") +
ggtitle("Relative difference in SD") +
theme(axis.text.x = element_text(angle=45, hjust=1))grpvitT <- as_labeller(c(grow = "Growth", mort="Mortality",
rec="Recruitment",
quadrat = "Space", `sp:time`="Species x Time",
`quadrat:time`="Space x Time", `quadrat:sp`= "Species x Space",
sp="Species", Residual = "Residual",
time= "Time"))
pm <- dif %>%
ggplot(aes(x=p.rare, y=VPC)) +
geom_point(aes(col=fplot, shape=group), size=3) +
facet_grid(vital~term, labeller = grpvitT) +
scale_x_continuous(breaks=c(50,55,60,65)) +
scale_color_discrete(labels=c("Barro Colorado", "Fushan",
"Lambir", "Luquillo", "Pasoh"))+
ylab("Difference in VPC to original data") +
xlab("Proportion of rare species (%)") +
geom_hline(yintercept = 0, linetype= "dotted") +
theme(legend.position="bottom",
panel.spacing.y = unit(0.8, "cm"),
panel.background = element_rect(fill="white", color="black"))+
labs(tag="b)")
pm# png("figs/FIG_S4.3b_exclude_regroup_rare_VPC_time.png", height = 600, width=1010, res=100)
# pm
# dev.off()pmr <- difr %>%
ggplot(aes(x=p.rare, y=stdev)) +
geom_point(aes(col=fplot, shape=group), size=3) +
facet_grid(vital~term, labeller = grpvitT) +
scale_x_continuous(breaks=c(64,68,72)) +
scale_color_discrete(labels=c("Barro Colorado", "Fushan",
"Lambir", "Luquillo", "Pasoh"))+
ylab("Proportional difference in SD to original data") +
xlab("Proportion of rare species (%)") +
geom_hline(yintercept = 0, linetype= "dotted") +
theme(legend.position="bottom",
panel.spacing.y = unit(0.8, "cm"),
panel.background = element_rect(colour="lightgray"))
pmr# png("figs/exclude_regroup_rare_SD_time_relativeDIFF.png", height = 600, width=1010, res=100)
# pmr
# dev.off()